Biomolecular interactions determine how transcription factors recognize their DNA binding sites, how proteins interact with each other, and consequently how a biological system functions. Since many biological molecules bear considerable electric charge, electrostatic interactions are among the most important when studying biomolecular interactions. However, electrostatic interactions in biological systems are difficult to calculate accurately in practice. Aside from the significant charges carried by biomolecules such as DNA and proteins, the solvent itself, namely, water produces considerable electrostatic effects. Furthermore, hydrogen bonds, known to be involved in helix formation in both DNA and proteins, are essentially electrostatic in origin. Indeed, it seems that electrostatic effects often drive the physical-chemical processes in biological systems and, thereby, determine biological functions. Therefore, any attempt to perform molecular dynamics (MD) simulations of biological systems will require an adequate description of these electrostatic forces. Previously, we have developed a rigorous surface charge method (SCM) to calculate the crucial electrostatic forces in a biomolecular system. Exact analytical results have been obtained for systems with sufficient symmetry. For example we developed methods to compute the electrostatic forces for a biomolecular system in which the atoms are represented by spheres. The method is rigorous in the context of the model and the accuracy can be tuned to any desired level. Our efforts in the past years have been to implement this exact formalism in a numerical code aiding realistic inclusion of polarization effect in calculating electrostatic interactions among biomolecules. The numerical code, although founded on an exact formalism, is slow due to the computation of many Wigner rotation matrix elements. This year, we have made progress in this direction. We are able to express the interaction energy among dielectric spheres in a closed form bypassing the use Wigner rotation matrices. This new but equivalent formalism was tested for systems of up to six dielectric spheres. It was found that the speed is increased by 20 fold or more. This important result makes possible molecular dynamics simulations with precise electrostatic energy. It is recently published in Physical Review E.